! Project TWS Matrixexponentiële
!
!
! MODULE IO
!
! Deze module bevat routines en functies om volle en circulante matrices weg te schrijven en in te lezen.
!
! Het formaat voor volle matrices is:
! GROOTTE
! RIJ1
! RIJ2
! ...
!
! Het formaat voor circulante matrices is:
! cGROOTTE
! RIJ1(1)
! RIJ2(1)
! ...
!
! @author Roel Matthysen
! @version 2.0

module IO

    implicit none
    private

    integer, parameter :: doubleprecision = kind(0.d0)    
    public :: printmatrix, readmatrix, iscirculant, readmatrixc, displayhelp, printresult
    
    
    ! PRINTMATRIX
    !
    ! Deze interface zorgt ervoor dat alle printmethodes met dezelfde naam kunnen worden opgeroepen.
    ! Deze printmethodes drukken een reële/complexe matrix/vector af op de errorstream.
    ! Deze methodes zijn bedoeld voor debugging.
    !
    ! @since    2.0
    interface printmatrix
        module procedure printmatrixr
        module procedure printmatrixc
        module procedure printvectorr
        module procedure printvectorc
    end interface
    
    ! PRINTRESULT
    !
    ! Deze interface zorgt ervoor dat printresult de formaten voor volle en circulante matrices volgt.
    !
    ! @since    2.0
    interface printresult
        module procedure printresultfull
        module procedure printresultcirc
    end interface
    
    contains
        
        ! Deze functie leest een volle matrix in van een stream, en geeft het resultaat terug in een variabele
        !
        ! @param     n                    De grootte van de in te lezen matrix.
        ! 
        ! @param     optunit (optional)   Dit argument bepaalt de stream waarvan de matrix wordt ingelezen.
        !                                 Indien dit argument achterwege wordt gelaten wordt de standaardinvoer     
        !                                 gebruikt
        ! @param     advanceone(optional) Dit argument bepaalt of de eerste regel moet worden overgeslagen of niet.
        !                                 Indien dit argument achterwege wordt gelaten wordt vanaf de eerste regel ingelezen.
        ! @result    A                    Een tweedimensionale doubleprecision matrix
        !
        ! @since     1.0
        function readmatrix(n, optunit,advanceone) result(A)
            real(doubleprecision),allocatable :: A(:,:)
            integer :: i,j,n,inputunit
            integer, optional :: optunit
            logical, optional :: advanceone
            allocate(A(n,n))
            if(present(optunit)) then
                inputunit = optunit
            else
                inputunit = 5
            endif
            if(present(advanceone)) then
                if(advanceone) then
                    read(inputunit,*)
                endif
            endif
            read(unit=inputunit,fmt=*) ((A(i,j),j=1,n),i=1,n)
        end function
        
        ! Deze functie leest een circulante matrix in van een stream, en geeft het resultaat 
        ! terug in een variabele
        !
        ! @param     n                    De grootte van de in te lezen matrix.
        !
        ! @param     optunit (optional)   Dit argument bepaalt de stream waarvan de matrix wordt ingelezen.
        !                                 Indien dit argument achterwege wordt geleten wordt de standaardinvoer     
        !                                 gebruikt    
        ! @param     advanceone(optional) Dit argument bepaalt of de eerste regel moet worden overgeslagen of niet.
        !                                 Indien dit argument achterwege wordt gelaten wordt vanaf de eerste regel ingelezen.
        ! @result    A                    Een ééndimensionale doubleprecision matrix
        !
        ! @since     2.0
        function readmatrixc(n, optunit, advanceone) result(A)
            real(doubleprecision),allocatable :: A(:)
            integer :: i,n,inputunit
            logical,optional :: advanceone
            integer, optional :: optunit
            allocate(A(n))
            if(present(optunit)) then
                inputunit = optunit
            else
                inputunit = 5
            endif
            if(present(advanceone)) then
                if(advanceone) then
                    read(inputunit,*)
                endif
            endif
            do i=1,n
                read(inputunit,*) A(i)
            enddo
        end function
    
    
        ! Deze subroutine bepaalt of de matrix die ingelezen wordt van de stream circulant is. De 
        ! 
        ! @param     inputunit            Dit argument bepaalt de stream waarvan de matrix wordt ingelezen.
        !                                    
        ! @result    circulant            Een logical die aanduidt of de matrix al dan niet circulant is.
        !
        ! @result    n                    De dimensie van de matrix (voor zowel het volle als het circulante geval).
        !
        ! @since     2.0
        subroutine iscirculant(inputunit, circulant, n)
            logical, intent(out) :: circulant
            integer, intent(in) :: inputunit
            integer :: n
            character (len=32) :: header
            read(inputunit,'(a)') header
            header = adjustl(header)
            circulant = (header(1:1)=='c')
            if(circulant) then
                header = header(2:32)
            endif
            header=trim(header);
            read(header,'(i32)') n
        end subroutine
        
        
        ! PRINTMATRIXR
        !
        ! Deze routine print een reële tweedimensionele matrix naar een stream met eventueel een header.
        !
        ! @param    A                    De te printen tweedimensionale reële matrix
        !
        ! @param    header(optional)     De toe te voegen header boven de matrix. Indien niet aanwezig wordt er niets
        !                                geprint.
        !
        ! @param    optunit(optional)    De stream waarop geprint moet worden. Indien niet aanwezig wordt op de 
        !                                standaard error stream (0) geschreven.
        ! @since    2.0
        subroutine printmatrixr(A, header, optunit)
            character(len=*),optional :: header
            integer, optional :: optunit
            real(doubleprecision), dimension(:,:) :: A
            integer :: i, outputunit = 0
            if(present(optunit)) outputunit = optunit
            if(present(header)) then
                write(outputunit,*) header
            endif
            do i=1,size(A,1)
                write(outputunit,*) A(i,:)
            end do
        end subroutine
    
        
        ! PRINTMATRIXR
        !
        ! Deze routine print een complexe tweedimensionele matrix naar een stream met eventueel een header.
        !
        ! @param    A                    De te printen tweedimensionale complexe matrix
        !
        ! @param    header(optional)     De toe te voegen header boven de matrix. Indien niet aanwezig wordt er niets
        !                                geprint.
        !
        ! @param    optunit(optional)    De stream waarop geprint moet worden. Indien niet aanwezig wordt op de 
        !                                standaard error stream (0) geschreven.
        ! @since    2.0
        subroutine printmatrixc(A, header, optunit)
            character(len=*), optional :: header
            integer, optional :: optunit
            complex(doubleprecision), dimension(:,:) :: A
            integer :: i, outputunit = 0
            if(present(optunit)) outputunit = optunit
            if(present(header)) then
                write(outputunit,*) header
            endif
            do i=1,size(A,1)
                write(outputunit,*) A(i,:)
            end do
        end subroutine
    
        ! PRINTVECTORR
        !
        ! Deze routine print een reële ééndimensionele matrix naar een stream met eventueel een header.
        !
        ! @param    A                    De te printen ééndimensionale reële matrix
        !
        ! @param    header(optional)     De toe te voegen header boven de matrix. Indien niet aanwezig wordt er niets
        !                                geprint.
        !
        ! @param    optunit(optional)    De stream waarop geprint moet worden. Indien niet aanwezig wordt op de 
        !                                standaard error stream (0) geschreven.
        ! @since    2.0
        subroutine printvectorr(A, header, optunit)
            character(len=*), optional :: header
            integer, optional :: optunit
            real(doubleprecision), dimension(:) :: A
            integer :: i, outputunit = 0
            if(present(optunit)) outputunit = optunit
            if(present(header)) then
                write(outputunit,*) header    
            endif
            do i=1,size(A,1)
                write(outputunit,*) A(i)
            end do
        end subroutine
        
        
        ! PRINTVECTORR
        !
        ! Deze routine print een complexe ééndimensionele matrix naar een stream met eventueel een header.
        !
        ! @param    A                    De te printen ééndimensionale complexe matrix
        !
        ! @param    header(optional)     De toe te voegen header boven de matrix. Indien niet aanwezig wordt er niets
        !                                geprint.
        !
        ! @param    optunit(optional)    De stream waarop geprint moet worden. Indien niet aanwezig wordt op de 
        !                                standaard error stream (0) geschreven.
        ! @since    2.0
        subroutine printvectorc(A, header, optunit)
            character(len=*), optional :: header
            integer, optional :: optunit
            complex(doubleprecision), dimension(:) :: A
            integer :: i, outputunit = 0
            if(present(optunit)) outputunit = optunit
            if(present(header)) then
                write(outputunit,*) header
            endif    
            do i=1,size(A,1)
                write(outputunit,*) A(i)
            end do
        end subroutine
        
        ! PRINTRESULTFULL
        !
        ! Deze routine print een volle matrix naar een stream volgens het formaat voor volle matrices.
        !
        ! @param    A                    De te printen volle matrix (reëel).
        !
        ! @param    optunit(optional)    De stream waarop geprint moet worden. Indien niet aanwezig wordt op de 
        !                                standaard uitvoer stream (6) geschreven.
        ! @since    2.0
        subroutine printresultfull(A, optunit)
            real(doubleprecision), dimension(:,:) :: A
            character(64) :: header
            integer, optional :: optunit
            integer :: outputunit = 6
            if(present(optunit)) outputunit = optunit
            write(header,*) size(A,1)
            call printmatrixr(A,header,outputunit)
        end subroutine
        
        ! PRINTRESULTFULL
        !
        ! Deze routine print een circulante matrix naar een stream volgens het formaat voor circulante matrices.
        !
        ! @param    A                    De te printen circulante matrix (reëel en voorgesteld door de eerste kolom).
        !
        ! @param    optunit(optional)    De stream waarop geprint moet worden. Indien niet aanwezig wordt op de 
        !                                standaard uitvoer stream (6) geschreven.
        ! @since    2.0
        subroutine printresultcirc(A, optunit)
            real(doubleprecision), dimension(:) :: A
            character(64) :: header
            integer, optional :: optunit
            integer :: outputunit = 6
            if(present(optunit)) outputunit = optunit
            write(header,'(i64)') size(A,1)
            header = 'c' // adjustl(header)
            call printvectorr(A,header,outputunit)
        end subroutine
        
        ! DISPLAYHELP
        !
        ! Deze routine geeft de help-boodschap weer op de standard error stream.
        subroutine displayhelp()
            write(0,'(a)') ' '
            write(0,'(a)') 'HELP MATRIXEXP'
            write(0,'(a)') ' '
            write(0,'(a)') '    Dit programma rekent voor een gegeven matrix de matrixexponentiële'             
            write(0,'(a)') '    uit. Indien geen invoer en uitvoer meegegeven worden wordt er'                 
            write(0,'(a)') '    ingelezen en uitgeschreven op de standaardstream. De invoer moet'                 
            write(0,'(a)') '    op de eerste lijn de dimensie van de matrix hebben, en op de'             
            write(0,'(a)') '    volgende lijnen de data. Het programma werkt voor zowel volle'
            write(0,'(a)') '    als circulante matrices, indien de input juist geformatteerd is'
            write(0,'(a)') ' '
            write(0,'(a)') 'OPTIONS'
            write(0,'(a)') ' '
            write(0,'(a)') '    --version           Toont informatie over versie en auteur'
            write(0,'(a)') ' '
            write(0,'(a)') '    --help              Toont deze boodschap'
            write(0,'(a)') ' '
            write(0,'(a)') '    --algorithm [ALG]   Keuze voor het algoritme. De keuzes voor ALG zijn:'
            write(0,'(a)') '                        taylor:  Op basis van reeksontwikkeling'
            write(0,'(a)') '                        eig:        Op basis van de eigenwaardenontbinding'
            write(0,'(a)') '                        eigqr:      Op basis van de eigenwaardenontbinding'
            write(0,'(a)') '                                    en QR-ontbinding'
            write(0,'(a)') ' '
            write(0,'(a)') '    --verbose           Zet het printen van debugging-informatie aan.'
        end subroutine
end module
